Chapter4.R

library(tidycensus)

ga_wide <- get_acs(
  geography = "county",
  state = "Georgia",
  variables = c(medinc = "B19013_001",
                medage = "B01002_001"),
  output = "wide",
  year = 2020
)

## Getting data from the 2016-2020 5-year ACS

## Using FIPS code '13' for state 'Georgia'

library(tidyverse)
options(scipen = 999)

ggplot(ga_wide, aes(x = medincE)) + 
  geom_histogram()

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

ggplot(ga_wide, aes(x = medincE)) + 
  geom_histogram(bins = 15)

ggplot(ga_wide, aes(y = medincE)) + 
  geom_boxplot()

ggplot(ga_wide, aes(x = medageE, y = medincE)) + 
  geom_point()

ggplot(ga_wide, aes(x = medageE, y = medincE)) + 
  geom_point() + 
  geom_smooth(method = "lm")

## `geom_smooth()` using formula = 'y ~ x'

library(tidycensus)
library(tidyverse)

metros <-  get_acs(
  geography = "cbsa",
  variables = "DP03_0021P",
  summary_var = "B01003_001",
  survey = "acs1",
  year = 2019
) %>%
  slice_max(summary_est, n = 20)

## Getting data from the 2019 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using the ACS Data Profile

ggplot(metros, aes(x = NAME, y = estimate)) + 
  geom_col()

metros %>%
  mutate(NAME = str_remove(NAME, "-.*$")) %>%
  mutate(NAME = str_remove(NAME, ",.*$")) %>%
  ggplot(aes(y = reorder(NAME, estimate), x = estimate)) + 
  geom_col()

metros %>%
  mutate(NAME = str_remove(NAME, "-.*$")) %>%
  mutate(NAME = str_remove(NAME, ",.*$")) %>%
  ggplot(aes(y = reorder(NAME, estimate), x = estimate)) + 
  geom_col() +  
  theme_minimal() + 
  labs(title = "Public transit commute share", 
       subtitle = "2019 1-year ACS estimates", 
       y = "", 
       x = "ACS estimate", 
       caption = "Source: ACS Data Profile variable DP03_0021P via the tidycensus R package") 

library(scales)

metros %>%
  mutate(NAME = str_remove(NAME, "-.*$")) %>%
  mutate(NAME = str_remove(NAME, ",.*$")) %>%
  ggplot(aes(y = reorder(NAME, estimate), x = estimate)) + 
  geom_col(color = "navy", fill = "navy", 
           alpha = 0.5, width = 0.85) +  
  theme_minimal(base_size = 12, base_family = "Verdana") + 
  scale_x_continuous(labels = label_percent(scale = 1)) + 
  labs(title = "Public transit commute share", 
       subtitle = "2019 1-year ACS estimates", 
       y = "", 
       x = "ACS estimate", 
       caption = "Source: ACS Data Profile variable DP03_0021P via the tidycensus R package") 

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

ggsave("metro_transit.png")

## Saving 7 x 5 in image

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

ggsave(
  filename = "metro_transit.png",
  path = "~/images",
  width = 8,
  height = 5,
  units = "in",
  dpi = 300
)

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

maine <- get_decennial(
  state = "Maine",
  geography = "county",
  variables = c(totalpop = "P1_001N"),
  year = 2020
) %>%
  arrange(desc(value))

## Getting data from the 2020 decennial Census

## Using FIPS code '23' for state 'Maine'

## Using the PL 94-171 Redistricting Data summary file

maine_income <- get_acs(
  state = "Maine",
  geography = "county",
  variables = c(hhincome = "B19013_001"),
  year = 2020
) %>%
  mutate(NAME = str_remove(NAME, " County, Maine"))

## Getting data from the 2016-2020 5-year ACS

## Using FIPS code '23' for state 'Maine'

ggplot(maine_income, aes(x = estimate, y = reorder(NAME, estimate))) + 
  geom_point(size = 3, color = "darkgreen") + 
  labs(title = "Median household income", 
       subtitle = "Counties in Maine", 
       x = "", 
       y = "ACS estimate") + 
  theme_minimal(base_size = 12.5) + 
  scale_x_continuous(labels = label_dollar())

maine_income %>% 
  arrange(desc(moe))

## # A tibble: 16 × 5
##    GEOID NAME         variable estimate   moe
##    <chr> <chr>        <chr>       <dbl> <dbl>
##  1 23023 Sagadahoc    hhincome    68039  4616
##  2 23015 Lincoln      hhincome    58125  3974
##  3 23027 Waldo        hhincome    58034  3482
##  4 23007 Franklin     hhincome    51630  2948
##  5 23021 Piscataquis  hhincome    42083  2883
##  6 23025 Somerset     hhincome    45382  2694
##  7 23009 Hancock      hhincome    58345  2593
##  8 23013 Knox         hhincome    57794  2528
##  9 23017 Oxford       hhincome    49761  2380
## 10 23003 Aroostook    hhincome    43791  2306
## 11 23029 Washington   hhincome    44847  2292
## 12 23031 York         hhincome    68932  2239
## 13 23011 Kennebec     hhincome    55368  2112
## 14 23001 Androscoggin hhincome    55002  2003
## 15 23019 Penobscot    hhincome    52128  1836
## 16 23005 Cumberland   hhincome    76014  1563

ggplot(maine_income, aes(x = estimate, y = reorder(NAME, estimate))) + 
  geom_errorbarh(aes(xmin = estimate - moe, xmax = estimate + moe)) + 
  geom_point(size = 3, color = "darkgreen") + 
  theme_minimal(base_size = 12.5) + 
  labs(title = "Median household income", 
       subtitle = "Counties in Maine", 
       x = "2016-2020 ACS estimate", 
       y = "") + 
  scale_x_continuous(labels = label_dollar())

years <- 2005:2019
names(years) <- years

deschutes_value <- map_dfr(years, ~{
  get_acs(
    geography = "county",
    variables = "B25077_001",
    state = "OR",
    county = "Deschutes",
    year = .x,
    survey = "acs1"
  )
}, .id = "year")

## Getting data from the 2005 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2006 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2007 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2008 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2009 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2010 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2011 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2012 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2013 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2014 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2015 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2016 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2017 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2018 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

## Getting data from the 2019 1-year ACS

## The 1-year ACS provides data for geographies with populations of 65,000 and greater.

## Using FIPS code '41' for state 'OR'

## Using FIPS code '017' for 'Deschutes County'

ggplot(deschutes_value, aes(x = year, y = estimate, group = 1)) + 
  geom_line() + 
  geom_point()

ggplot(deschutes_value, aes(x = year, y = estimate, group = 1)) + 
  geom_ribbon(aes(ymax = estimate + moe, ymin = estimate - moe), 
              fill = "navy",
              alpha = 0.4) + 
  geom_line(color = "navy") + 
  geom_point(color = "navy", size = 2) + 
  theme_minimal(base_size = 12) + 
  scale_y_continuous(labels = label_dollar(scale = .001, suffix = "k")) + 
  labs(title = "Median home value in Deschutes County, OR",
       x = "Year",
       y = "ACS estimate",
       caption = "Shaded area represents margin of error around the ACS estimate")

utah <- get_estimates(
  geography = "state",
  state = "UT",
  product = "characteristics",
  breakdown = c("SEX", "AGEGROUP"),
  breakdown_labels = TRUE,
  year = 2019
) 

## Using FIPS code '49' for state 'UT'

utah_filtered <- filter(utah, str_detect(AGEGROUP, "^Age"), 
                        SEX != "Both sexes") %>%
  mutate(value = ifelse(SEX == "Male", -value, value))

ggplot(utah_filtered, aes(x = value, y = AGEGROUP, fill = SEX)) + 
  geom_col()

utah_pyramid <- ggplot(utah_filtered, 
                       aes(x = value, 
                           y = AGEGROUP, 
                           fill = SEX)) + 
  geom_col(width = 0.95, alpha = 0.75) + 
  theme_minimal(base_family = "Verdana", 
                base_size = 12) + 
  scale_x_continuous(
    labels = ~ number_format(scale = .001, suffix = "k")(abs(.x)),
    limits = 140000 * c(-1,1)
  ) + 
  scale_y_discrete(labels = ~ str_remove_all(.x, "Age\\s|\\syears")) + 
  scale_fill_manual(values = c("darkred", "navy")) + 
  labs(x = "", 
       y = "2019 Census Bureau population estimate", 
       title = "Population structure in Utah", 
       fill = "", 
       caption = "Data source: US Census Bureau population estimates & tidycensus R package")

utah_pyramid

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call.graphics(C_text, as.graphicsAnnot(x$label), x$x, x$y, : font
## family not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

## Warning in grid.Call(C_textBounds, as.graphicsAnnot(x$label), x$x, x$y, : font family
## not found in Windows font database

housing_val <- get_acs(
  geography = "tract", 
  variables = "B25077_001", 
  state = "OR", 
  county = c(
    "Multnomah", 
    "Clackamas", 
    "Washington",
    "Yamhill", 
    "Marion", 
    "Columbia"
  ),
  year = 2020
)

## Getting data from the 2016-2020 5-year ACS

## Using FIPS code '41' for state 'OR'

## Using FIPS code '051' for 'Multnomah County'

## Using FIPS code '005' for 'Clackamas County'

## Using FIPS code '067' for 'Washington County'

## Using FIPS code '071' for 'Yamhill County'

## Using FIPS code '047' for 'Marion County'

## Using FIPS code '009' for 'Columbia County'

housing_val2 <- separate(
  housing_val, 
  NAME, 
  into = c("tract", "county", "state"), 
  sep = ", "
)

housing_val2 %>%
  group_by(county) %>%
  summarize(min = min(estimate, na.rm = TRUE), 
            mean = mean(estimate, na.rm = TRUE), 
            median = median(estimate, na.rm = TRUE), 
            max = max(estimate, na.rm = TRUE))

## # A tibble: 6 × 5
##   county               min    mean median     max
##   <chr>              <dbl>   <dbl>  <dbl>   <dbl>
## 1 Clackamas County   62800 449941. 426700  909000
## 2 Columbia County   218100 277591. 275900  362200
## 3 Marion County      48700 270969. 261200  483500
## 4 Multnomah County  192900 455706. 425950 1033500
## 5 Washington County 221900 419618. 406100  769700
## 6 Yamhill County    230000 333816. 291100  545500

ggplot(housing_val2, aes(x = estimate)) + 
  geom_density()

## Warning: Removed 9 rows containing non-finite values (`stat_density()`).

ggplot(housing_val2, aes(x = estimate, fill = county)) + 
  geom_density(alpha = 0.3)

## Warning: Removed 9 rows containing non-finite values (`stat_density()`).

ggplot(housing_val2, aes(x = estimate)) +
  geom_density(fill = "darkgreen", color = "darkgreen", alpha = 0.5) + 
  facet_wrap(~county) + 
  scale_x_continuous(labels = dollar_format(scale = 0.000001, 
                                            suffix = "m")) + 
  theme_minimal(base_size = 14) + 
  theme(axis.text.y = element_blank(), 
        axis.text.x = element_text(angle = 45)) + 
  labs(x = "ACS estimate",
       y = "",
       title = "Median home values by Census tract, 2015-2019 ACS")

## Warning: Removed 9 rows containing non-finite values (`stat_density()`).

library(ggridges)

ggplot(housing_val2, aes(x = estimate, y = county)) + 
  geom_density_ridges() + 
  theme_ridges() + 
  labs(x = "Median home value: 2016-2020 ACS estimate", 
       y = "") + 
  scale_x_continuous(labels = label_dollar(scale = .000001, suffix = "m"),
                     breaks = c(0, 500000, 1000000)) + 
  theme(axis.text.x = element_text(angle = 45))

## Picking joint bandwidth of 36200

## Warning: Removed 9 rows containing non-finite values (`stat_density_ridges()`).

library(ggbeeswarm)

ny_race_income <- get_acs(
  geography = "tract", 
  state = "NY",  
  county = c("New York", "Bronx", "Queens", "Richmond", "Kings"),
  variables = c(White = "B03002_003", 
                Black = "B03002_004", 
                Asian = "B03002_006",
                Hispanic = "B03002_012"), 
  summary_var = "B19013_001",
  year = 2020
) %>%
  group_by(GEOID) %>%
  filter(estimate == max(estimate, na.rm = TRUE)) %>%
  ungroup() %>%
  filter(estimate != 0)

## Getting data from the 2016-2020 5-year ACS

## Using FIPS code '36' for state 'NY'

## Using FIPS code '061' for 'New York County'

## Using FIPS code '005' for 'Bronx County'

## Using FIPS code '081' for 'Queens County'

## Using FIPS code '085' for 'Richmond County'

## Using FIPS code '047' for 'Kings County'

ggplot(ny_race_income, aes(x = variable, y = summary_est, color = summary_est)) +
  geom_quasirandom(alpha = 0.5) + 
  coord_flip() + 
  theme_minimal(base_size = 13) + 
  scale_color_viridis_c(guide = "none") + 
  scale_y_continuous(labels = label_dollar()) + 
  labs(x = "Largest group in Census tract", 
       y = "Median household income", 
       title = "Household income distribution by largest racial/ethnic group", 
       subtitle = "Census tracts, New York City", 
       caption = "Data source: 2016-2020 ACS")

## Warning: Removed 35 rows containing missing values.

library(geofacet)

us_pyramid_data <- get_estimates(
  geography = "state",
  product = "characteristics",
  breakdown = c("SEX", "AGEGROUP"),
  breakdown_labels = TRUE,
  year = 2019
) %>%
  filter(str_detect(AGEGROUP, "^Age"),
         SEX != "Both sexes") %>%
  group_by(NAME) %>%
  mutate(prop = value / sum(value, na.rm = TRUE)) %>%
  ungroup() %>%
  mutate(prop = ifelse(SEX == "Male", -prop, prop))

ggplot(us_pyramid_data, aes(x = prop, y = AGEGROUP, fill = SEX)) + 
  geom_col(width = 1) + 
  theme_minimal() + 
  scale_fill_manual(values = c("darkred", "navy")) + 
  facet_geo(~NAME, grid = "us_state_with_DC_PR_grid2",
            label = "code") + 
  theme(axis.text = element_blank(),
        strip.text.x = element_text(size = 8)) + 
  labs(x = "", 
       y = "", 
       title = "Population structure by age and sex", 
       fill = "", 
       caption = "Data source: US Census Bureau population estimates & tidycensus R package")

library(plotly)

ggplotly(utah_pyramid)